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Abstract 

We utilize classical molecular dynamics to study surface effects 
on the piezoelectric properties of ZnO nanowires as calculated under 
uniaxial loading. An important point to our work is that we have uti- 
lized two types of surface treatments, those of charge compensation 
and surface passivation, to eliminate the polarization divergence that 
otherwise occurs due to the polar (0001) surfaces of ZnO. In doing so, 
we find that if appropriate surface treatments are utilized, the elastic 
modulus and the piezoelectric properties for ZnO nanowires having a 
variety of axial and surface orientations are all reduced as compared to 
the bulk value as a result of polarization- reduction in the polar [0001] 
direction. The reduction in effective piezoelectric constant is found 
to be independent of the expansion or contraction of the polar (0001) 
surface in response to surface stresses. Instead, the surface polariza- 
tion and thus effective piezoelectric constant is substantially reduced 
due to a reduction in the bond length of the Zn-0 dimer closest to the 
polar (0001) surface. Furthermore, depending on the nanowire axial 
orientation, we find in the absence of surface treatment that the piezo- 
electric properties of ZnO are either effectively lost due to unphysical 
transformations from the wurtzite to non-piezoelectric d-BCT phases, 
or also become smaller with decreasing nanowire size. The overall 
implication of this study is that if enhancement of the piezoelectric 
properties of ZnO is desired, then continued miniaturization of square 
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or nearly square cross section ZnO wires to the nanometer scale is not 
likely to achieve this result. 

1 Introduction 

Piezoelectricity has long been a property of interest for bulk materials as 
it enables the direct conversion of mechanical strain into harvestable elec- 
trical energy |49l |6] . While the interest in bulk piezoelectric materials has 
existed for some time, there has recently been significant interest in study- 
ing the piezoelectric behavior and properties of nanomaterials [60] . Much of 
the interest has centered around ZnO, which was recently utilized by Wang 
et al ^ to generate electrical energy through application of bending 
deformation via an atomic force microscope (AFM). ZnO has proven to be 
a versatile choice for nanoscale piezoelectrics as it exhibits both semicon- 
ducting and piezoelectric properties [68j, because it can be fabricated in a 
wide range of nanometer shapes and geometries [67], and because it has 
the largest piezoelectric response of any tetrahedrally bonded semiconduc- 
tor [13]. Since the initial discovery in 2006, there have since emerged, though 
not without controversy [5] , a wide range of interesting applications involving 
ZnO [Ml in], GaN [55], and other nanostructures [IH El EH [59l [56] . 

In addition to the wide range of potential applications, recent experimen- 
tal [80j and computational [TQl IHl [161 132] work has suggested that due to 
nanoscale surface eflFects, ZnO nanostructures may exhibit different piezo- 
electric properties than bulk ZnO. These non-bulk piezoelectric properties 
may couple with the recent finding that ZnO nanostructures exhibit me- 
chanical properties, and specifically Young's modulus that also shows a clear 
size-dependence due to surface eflFects fT2\ l2l 148] to potentially enable ZnO 
nanowires (NWs) to produce more mechanical strain energy that can be con- 
verted through the piezoelectric eflFect into harvestable electrical energy than 
bulk ZnO. 

However, a key issue that has not been resolved is how surface effects 
impact the piezoelectric properties of ZnO NWs. In other words, will mak- 
ing ZnO NWs smaller lead to enhanced piezoelectric properties? We note 
that the NW geometry has been studied for other materials, for example 
using molecular dynamics (MD) for BTO [ZZl [ZHl [79] , for GaN NWs using 
ab initio techniques ^Lj, and also using recently developed analytical theo- 
ries [521 [731 [46l|38l [37] . The piezoelectric properties of ZnO nanostructures. 
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though excluding surface effects, have also been studied primarily using ab 
initio calculations |28l I4j ; the surface piezoelectric properties of ZnO were 
recently studied by [T6], though the effects on one-dimensional NWs were 
not considered. A recent MD study did consider ZnO nanobelts [42j, though 
only for the [0001] orientation in which the transverse surfaces are not the 
polar (0001) surfaces and in which surface treatment, as described in the fol- 
lowing paragraph, were not considered. The one-dimensional NW geometry 
is critical to study and understand because it is most often utilized in appli- 
cation [66j, where the NWs are subject to axial flU [75[ [72] , bending [68^ 65j 
or shear deformations [39]. Furthermore, ZnO NWs can be synthesized with 
a variety of axial and surface orientations [80], and cross sectional geome- 
tries [1], which will impact the piezoelectric properties in different fashions. 
Therefore, a comprehensive understanding of how surface effects impact the 
piezoelectric properties of ZnO NWs, and how the piezoelectric properties 
of ZnO NWs vary with different surface and axial orientations of ZnO NWs 
remains unresolved. It is the purpose of this work to shed insight into these 
issues, by virtue of classical MD simulations. 

A related, and important issue this work addresses is the effect of the 
treatment of the polar ZnO (0001) surfaces on the piezoelectric properties. 
Specifically, as previously discussed by and subsequently by other re- 
searchers [471 ESI [29] , when there is a dipole moment in the repeat unit nor- 
mal to the surface of an ionic crystal, the electrostatic energy diverges, and 
the surface energy goes to infinity. Because of this, there are typically three 
techniques that are employed in atomistic simulations to eliminate this effect: 
charge compensation (CC) surface reconstruction (SR) and 

surface passivation (SP) or adsorption [54j. These stabilization techniques 
are utilized because such reconstructions and passivation have been observed 
experimentally [SH [33l [I9] , and they have also been widely used in first prin- 
ciples calculations fT4^, ^61]. In contrast, they have rarely been utilized in 
classical MD simulations [26l HE] to avoid the divergence of the electrostatic 
potential. While it is crucial to adopt one of these surface treatments for 
electrostatic stabilization, such treatments have not been utilized in previous 
MD studies of the size-dependent elastic properties of ZnO [3^, or ab ini- 
tio studies of the piezoelectric properties of other ZnO NWs [HJ [I] . We will 
demonstrate the issues that arise in the electromechanical properties of ZnO 
if no surface treatment is undertaken. 
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2 Methods 



We utilized classical MD to study the piezoelectric properties of ZnO NWs. 
Specifically, we used the open source GROMACS 4.0 molecular simulation 
code [24] while employing the Buckingham potential of [10] to model the 
various Zn-0 interactions. The Binks potential has been widely utilized to 
study the mechanical deformation of ZnO NWs [31]. However, until recent 
work by the authors for both bulk ZnO [15j, and subsequently for the surfaces 
of ZnO ^], the performance of the Binks potential for the piezoelectric 
properties of ZnO had not been investigated. Both works found the accuracy 
of the classical Binks potential to be comparable to benchmark ab initio 
calculation results \13\ . 

The lattice parameters we used for ZnO were ao=3.2709 A, Co=5.1386 A 
and ii=0.3882. For the electrostatic interactions, we utilized the approach 
of [20] , who improved on the original work of [69j by ensuring that the elec- 
trostatic potential and force smoothly truncate at the cut-off radius, which 
results in stability for MD simulations [21j. The approach of ^20j, which 
enables the convergent calculation of the electrostatic energies and forces us- 
ing a finite cut-off distance, is needed for the present simulations due to the 
fact that the standard Ewald method assumes an infinite, periodic crystal 
which is certainly not the case here due to the surface-dominated NW ge- 
ometries. The errors introduced by using the Ewald summation as compared 
to the Wolf technique for finite-sized NWs were recently quantified by [23j. 
The specific parameters for the Fennell method that we utilized for ZnO 
were a = 3 nm~^, rc=l nm; these parameters were previously found to give 
convergent results for the piezoelectric properties of ZnO [16j. 

We considered nearly square cross section ZnO NWs with cross sectional 
lengths ranging from 2 to 4 nm. We did not consider NWs with cross sectional 
sizes smaller than 2 nm because at these small sizes, a transformation into 
either a nonpiezoelectric d-BCT lattice structure [30l [50l |63l |64] or a shell 
structure [31j occurred as was previously predicted using MD simulations. 
The specific combinations of axial and surface orientations we considered are 
illustrated in Fig. [T} with the NW sizes summarized in Table [T} where the 
[2110], [0110] and [0001] directions are always chosen to be parallel to the x, 
y and z axes. No periodic boundary conditions were utilized in any direction, 
which implies that a truly finite-sized NW geometry subject to surface effects 
was considered in the present work, and that the NW sizes listed in Table [T] 
are the actual sizes used for the MD simulations. 
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(a) 




(b) 



Figure 1: (a) Four atom unit cell (in green rectangular box) for the wurtzite 
crystal structure; (b) Three ZnO NWs considered in this work. From top, 
the axial orientations are along the [0110], [0001] and [2110] directions. The 
z-direction is chosen to be along the [0001] direction for all NW orientations. 
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Table 1: NW dimensions for all three orientations considered in Fig. [!} 
where A^^^, A^^ and A^^ represent the number of unit cells in each direction. 
The aspect ratio for each NW is chosen to be 4:1. 
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Figure 2: Illustration of fixed (left end (blue) and right end (green)) and free 
(red and gray) atoms for axial loading. The NW shown has a [0001] axial 
direction. 
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Figure 3: Illustration of surface passivation via (left) adsorption of an H 
atom on O terminated (OOOl), and (right) OH adsorbed on the Zn terminated 
(0001) surface leading to a 2 x 1 pattern for both surfaces, (red: O, grey: 
Zn, white: H) 



We performed MD simulations of tensile axial deformation. For the ten- 
sile loading, both ends of the NW were first allowed to relax dynamically to 
a new equilibrium length in response to surface stresses by using a Berendsen 
thermostat |9j for up to 400 ps depending on the NW cross sectional size. 
After the new equilibrium length was found, two unit cells at each end of the 
NW, as illustrated in Fig. [2| were held fixed while the NW was equilibrated 
using a Nose-Hoover thermostat [25] for 20 ps. After these two initial equilib- 
rium steps, the ends of the NW were displaced axially at strain increments 
of 0.25% and held fixed while the NW was relaxed for 20 ps. After each 
strain increment, both ends of the NW were held fixed while the remainder 
of the NW was dynamically equilibrated for 100 ps using the Nose-Hoover 
thermostat at a temperature of 300K. The loading was increased until an 
axial strain of 20% in tension was reached. 

For the surface treatments to avoid the electrostatic divergence, we note 
that the surface atoms on the (0001) polar surfaces of ZnO have three nearest 
neighbors instead of four as does a bulk Zn or O atom. For the SP treatment, 
we added an H atom to the 0-terminated surface, and added an OH molecule 
to the Zn terminated surface in order to saturate dangling bonds, which, as 
illustrated in Fig. [3| results in a 2x1 passivation on both the O and Zn- 
terminated polar (0001) surfaces. The 2x1 passivation is used in our work 
for multiple reasons. First, as previously mentioned it has been observed 
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Table 2: Buckingham parameters for the Zn-0 interactions from [10.J, and 
also for the H-0 interactions taken from [17] for the surface passivation. 
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experimentally [34] and has been used in previous DFT calculations p^H 
l6T] . Furthermore, it is also the smallest passivation pattern that we can use 
on the NW polar surfaces due to the relatively small cross-sectional sized 
NWs we consider in this work. The potential parameters for both the Zn-0 
interactions as modeled using the potential of [lOj, as well as the parameters 
for the H-0 interactions needed for the surface passivation as taken from [17j 
are listed in Table [2} This passivation is also realistic as it is common for 
the environment to contain some water or humidity; the effects of water 
on the elastic properties of ZnO have also recently been investigated [76j. 
We note the likelihood that the piezoelectric properties of ZnO will depend 
on the specific passivation that is utilized in computation, or that occurs 
experimentally. 

For the surface treatment using CC, the methodology is much more 
straightforward. Because each Zn and O atom on a polar (0001) surface 
has 75% of the neighbors of the corresponding bulk atom (i.e. 3 instead 
of 4), the charge of the top layer of Zn and O atoms is reduced to 75% of 
the formal charge from ±2e to ±1.5e [471 HBJ- The CC surface treatment 
can also be physically justified as enforcing partial covalence of the surface 
atoms as compared to bulk atoms. We also note that the SP and CC surface 
treatments can be used for other polar crystals [71 E^j • 

The key value of interest we will report is the change in polarization for 
the NW as a function of the applied mechanical deformation. This parameter 
is key for design of nanogenerators as a larger polarization is directly related 
to a higher output voltage [SU [22l [56] , and thus more electrical energy gen- 
eration [271 ITU E^. For the axial loading, we accomplish this by calculating, 
for each state of strain, the polarization of each unit cell, then summing over 
the entire NW to calculate the total NW polarization, where the unit cell is 



8 



defined by the group of four atoms in the green box in Fig. [T](a), i.e.: 

N 

^cell 



• ^ ^cell ■ 
1=1 j=l 



where Pceii is the polarization for a single unit cell, P3 is the polarization 
of the NW in the polar [0001] direction, is the total number of unit cells 
in the system, q is the charge on each atom and z is the coordinate of each 
atom. We note that the effective piezoelectric constants are calculated using 
the invariant definition of [58j : 

_eff _^dPs ^eff _ 1 dPs ff _ dPs . . 



3 Numerical Results 
3.1 Mechanical Properties 

We first discuss the effect that the different surface treatments (CC vs. SP) 
have on the elastic properties of ZnO NWs, where the results are summa- 
rized in Fig. [4} The Young's modulus for each orientation was calculated by 
normalizing by the bulk Young's modulus for each orientation, which were 
taken to be 156.2 GPa for the [2TT0] and [0110] orientations and 119.7 GPa 
for the [0001] orientation [35^. The results are consistent for all three NW 
orientations considered: the Young's modulus is highest when no surface 
treatment is performed (original), followed by SP followed by a substantial 
reduction for the CC case. Furthermore, the modulus is observed to increase 
with decreasing size for all three orientations for the original case, which is 
consistent with previous MD simulation reports on ZnO NW elastic proper- 
ties [31j, whereas a size-dependent decrease in Young's modulus is observed 
for both the CC and SP cases. The Young's modulus is lowest for the CC 
case because of the reduction in formal charge for the surface atoms, which 
leads to reduced interaction energies, forces and thus stiffness for the surface 
atoms, and for the [0001] and [0110] orientations leads to a Young's modulus 
that is smaller than the bulk value for the smallest NW sizes we considered. 
The modulus for the SP case is similar to the no treatment case, but is typ- 
ically slightly smaller and is found to be larger than the bulk value for all 
sizes considered. 
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Figure 4: Bulk-normalized Young's modulus for the three NW orientations 
and geometries summarized in Table [l] for CC, SP, and original (untreated) 
surface treatments. 



While the Young's modulus trends and values in Fig. |4]may seem surpris- 
ing, the atomistic origin of these trends can be observed as shown in Figs. 
[5] and |6| where the atomic structure before any axial loading is applied is 
shown for both the [0110] and [2110] orientations. Specifically, it is shown 
that for both the [0110] and [2110] orientations, if no surface treatment is 
performed, as has been the case in previous calculations [45l HH, |43j, the 
WZ lattice structure is unstable and transforms to a d-BCT structure. In 
fact, this transformation occurs during the initial relaxation phase of the 
simulation for all NW sizes we have considered, and therefore the Young's 
moduli for the [0110] and [2110] orientations in Fig. [Z] correspond to that of 
the d-BCT, and not WZ structure. We note that regardless of the surface 
treatment, no initial transformation occurs for the [0001] NWs because the 
polar direction is along the axial direction for this case. 

In contrast, the SP and CC cases for both the [0110] and [2110] orienta- 
tions result in a stable WZ structure for the NW sizes we have considered, 
which demonstrates if the polarization divergence due to the polar (0001) 
surfaces is not treated, the WZ lattice structure is not stable and trans- 
forms to the d-BCT structure. The transformation to the d-BCT phase will 
have significant ramifications on the piezoelectric constants, as we will dis- 
cuss shortly. Before continuing to that discussion, we note that after yield. 
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Figure 5: (a): Snapshot of the atomic configuration at zero tensile strain for 
SP, original and CC surface types for the [0110] orientation showing that the 
SP and CC NWs keep the original WZ lattice structure, while the original 
NW has transformed to a d-BCT phase, (b) Comparison of the (left) WZ 
lattice structure to the (right) d-BCT lattice structure (taken from the blue 
rectangular box in (a)), where two unit cells are chosen and highlighted in 
blue for comparison to show the typical rectangular structure of the d-BCT 
phase. 
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Figure 6: (a): Snapshot of the atomic configuration at zero tensile strain for 
SP, original and CC surface types for the [2110] orientation showing that the 
SP and CC NWs keep the original WZ lattice structure, while the original 
NW has transformed to a d-BCT phase, (b) Comparison of the (left) WZ 
lattice structure to the (right) d-BCT lattice structure (taken from the blue 
rectangular box in (a)), where two unit cells are chosen and highlighted in 
blue for comparison to show the typical rectangular structure of the d-BCT 
phase. 
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the NWs transform into a non-piezoelectric structure [30l EOj [63l Ell 131] . 
Specifically, the [0110] NW transforms to a hexagonal [32] phase, while the 
[0001] orientation transforms to the d-BCT phase [3j. The transformation to 
a non-piezoelectric phase can also be observed by the post-yield behavior in 
the polarization vs. strain curves in Figs. [7[ [8]and|9} which we discuss next 
in further detail. 

3.2 Piezoelectric Constants 

The polarization vs. strain for all three NW orientations is shown in Figs. 
[7| [8] and [9] The first issue to point out is that, if CC or SP is utilized 
for the [0110] and [2110] orientations in Figs, [t] and [sj the polarization is 
linearly dependent on strain until yield, which occurs around 10% tensile 
strain for both orientations. However, as seen in Figs. [7|(b) and [8](b), if 
no surface treatment is utilized, the slope of the polarization vs. strain 
curve varies significantly even at very small amounts of applied tensile strain. 
Furthermore, due to the transformation from the WZ to non-piezoelectric d- 
BCT structure, as shown in Figs. [5](b) and[6](b), the polarization vs. strain 
curves are quite noisy and do not converge with decreasing size, which is a 
direct result of the non-piezoelectric d-BCT structure. We note that for the 
[2110] orientation, the relevant effective piezoelectric constant is 63^, while 
for the [0110] orientation, the relevant effective piezoelectric constant is e^. 

The [0001] orientation exhibits a different polarization vs. strain response 
than the [0110] and [2110] orientations because, as previously discussed, the 
WZ phase is stable without any surface treatment. Therefore, as shown in 
Fig. [9| the slope of the polarization vs. strain curve is linear even for the 
original surface case. However, it is clear that the slope for the original 
surface is significantly higher than the slopes for the CC and SP cases. A 
simple argument as to why the original, untreated surface leads (incorrectly) 
to a larger piezoelectric constant can be given as follows: the dipole moment 
for the original NW with untreated surfaces can be estimated as Moriginai = 
2 * 2e * rfi A^, where di and d2 are the layer distances that are related by di = 
d2 Q , N is number of unit cells in the polar [0001] direction. The dipole 
moment of the CC NWs is Mcc = -2*2e*rf2(A-l)+2*1.5e(Arfi + (A-l)rf2), 
so ^^w^ - ^ = 1-475. This explains why in Table ji] the original is 
larger than from CC and SP, and demonstrates the necessity of CC and 
SP to obtain the more reasonable effective piezoelectric constants. 
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Figure 7: Polarization vs. strain for axial loading along the [2110] direc- 
tion for different NW sizes. Inset (a) Size-dependent effective piezoelectric 
coefficient 63^; (b) Polarization vs. strain for original (untreated surface) 
NW. 
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Figure 8: Polarization vs. strain for axial loading along the [0110] direc- 
tion for different NW sizes. Inset (a) Size-dependent effective piezoelectric 
coefficient Cgf ; (b) Polarization vs. strain for original (untreated surface) 
NW. 
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Figure 9: Polarization vs strain for axial loading along the [0001] direction 
for different NW sizes. Inset (a) Zoom in to the small strain regime; (b) 
Size-dependent effective piezoelectric coefficient e^. 



A summary of the effective piezoelectric constants for all orientations 
and sizes is given in Table [3| where a comparison of the bulk piezoelectric 
constants as calculated for the [lOj potential are given for reference. As can 
be seen, there is a decrease in effective piezoelectric constant with decreasing 
size if the CC and SP are utilized. Our previous MD and DFT study [16j also 
found that the effective piezoelectric coefficients of ZnO thin film decrease as 
the film thickness decreases. The decrease in Cgf is less dramatic, with the 
reduction reaching 16.6% for the smallest NW sizes considered. In contrast, 
the reduction in and is more dramatic, reaching 80.0% for the smallest 
NW sizes considered for egf . As found before for the Young's modulus, the 
reduction in the piezoelectric constants is generally larger if CC is utilized as 
compared to SP. 

We now address the mechanism underlying the smaller piezoelectric con- 
stants that we have found for the CC and SP surface treatments. The ap- 
proach we take, similar to previous works [8l [73 [T^, is to analyze the po- 
larization on a unit-cell basis through the polar [0001] direction of the NWs. 
For the [2110] and [0110] orientated NWs, the variation in unit cell polar- 



ization through the NW polar [0001] direction is shown in Fig. 10, where 
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Table 3: Summary of effective piezoelectric constants (units of C/m^) for the 
different NW sizes and orientations from Table [T] Comparison with the bulk 
piezoelectric constants from [15j provided for reference. Labeling for NWs is 
the same as in Table [T] for consistency. 
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the polarization at the surfaces corresponds to the polar (0001) surfaces. As 
shown in Fig. [T0| the unit cell polarization at the surface is reduced by more 
than 10% for the CC surface treatment for both NW orientations, while the 
surface unit cell polarization reduction is smaller, i.e. less than 5% for all 
SP-oriented NW sizes. The polarization reduction is largest at the polar 
surface for both CC and SP, then converges to the bulk value as the interior 
of the NW is reached. 

The atomistic deformation leading to the reduction in polarization for 
the surface unit cells as shown in Fig. [T0| and thus the reduction in effective 
piezoelectric constant as shown in Table |3] is shown in Fig. [TT} which shows 
a snapshot of a surface unit cell for the [2110] oriented ZnO NW with size 
A3, with atomic displacements resulting from both CC and SP surface treat- 
ments. It is at first glance surprising that both surface treatments lead to 
decreases in the effective piezoelectric constants with decreasing NW cross 
sectional size because the surface contracts for CC and expands for SP in the 
[0001] direction in response to surface stresses as shown in Fig. 11, However, 
the important parameter for the polarization is not the absolute displacement 
of atoms near the surface, but instead the relative displacements between the 
Zn-0 atoms that comprise each of the two dimers in the surface unit cell, as 
can be seen through inspection of Eq. [TJ Specifically, the relative distance 
between the Zn-0 dimer closest to the surface in Fig. [TT]is -0.0107 nm and 
-0.0018 nm for the CC and SP surface treatments, respectively. Similarly, 
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the bond length change between the Zn-0 dimer one dimer into the bulk is 
0.0063 nm and -0.0009 nm for the CC and SP surface treatments. It is thus 
clear that while the surface atoms show different relaxations for the CC and 
SP surface treatments, in both cases there is a decrease in distance between 
the Zn-0 dimer at the surface which is significantly larger than the distance 
change for the Zn-0 dimer that is one dimer into the bulk, and furthermore 
the bond length decrease is much larger for CC than SP. This relative de- 
crease in surface dimer bond length explains the decrease in polarization in 
Fig. [lOj and thus effective piezoelectric constant in Table [3| and also why 
the decrease in effective piezoelectric constant in Table |3] is much greater for 
the CC than SP surface treatments. 

To further investigate the validity of our calculated piezoelectric constants 
for the CC and SP surface treatments, we compare against existing DFT re- 
sults for the effective piezoelectric constants of ZnO NWs. Specifically, DFT 
calculations by [TOj and [36j also found a decrease in effective piezoelectric 
constants and with decreasing NW size, which is the same trend as 
found in the present work. We note that the comparison is not exact, as the 
NWs in the DFT calculations had a hexagonal cross section oriented in the 
[0001] direction that was enclosed by (0110) surfaces. However, [36j studied 
NW diameters from 3.1 to 0.4 nm and reported a decrease in from 1.5 
to 1.31 C/m^, for a reduction of 14.5%. In the present work, the reduction 
in the effective piezoelectric constant of our original, CC and SP rectan- 
gular NWs with cross sectional sizes from 2-4 nm oriented along the [0001] 
direction were found to be 14.9%, 9.13% and 13.0%, respectively. 

A final, but very important question to address is whether the CC and 
SP surface treatments will necessarily lead to a decrease in the effective 
piezoelectric constants of the NWs due to the fact that they reduce charge, 
and thus polarization at the NW surfaces. While our results did in fact 
show a reduction in piezoelectric constant with decreasing NW size for all 
NW geometries and orientations considered, other literature results suggest 
that this need not be the case. Specifically, we note the recent work of [Ij, 
who studied the piezoelectric properties of GaN and ZnO NWs, albeit with 
a hexagonal cross section as compared to the nearly square cross sections 
considered in the present work. They also observed a charge and polarization 
reduction with decreasing NW size due to surface effects. However, because of 
the strong contraction of the transverse surfaces, the reduction in volume (see 
Eq. ([T])), of the NW becomes more important than the reduction in surface 
charge, leading to a predicted increase in the effective piezoelectric constant 
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^ layer* ^ 



(a) 




Figure 10: Varia tion in uni t cell polarization through the NW [0001] direction 
for both the (a) [2110] and (b) [0110] axial orientations, which demonstrates 
the reduction in surface unit cell polarization as compared to the bulk. The 
polarization of each unit cell is normalized by the polarization of a bulk unit 
ceU. 
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Figure 11: Side view of the surface unit cell (four atoms, or two Zn-0 dimers) 
for the [2110] oriented ZnO NW with size A3. The red and blue arrows 
show the atomic displacements for CC and SP, respectively, with the values 
labeled in the corresponding boxes. The image shows that while the surface 
atoms contract for CC, but expand for SP in the [0001] direction, the relative 
distance between the Zn and O atoms that comprise the dimer nearest to 
the surface becomes smaller for both CC and SP. 
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for very small (< 2 nm diameter) NWs. These results also suggest that the 
cross sectional geometry may play a critical role in determining the size- 
dependence of the piezoelectric constant for NWs, as our preliminary studies 
also show an increase in effective piezoelectric constant with decreasing size 
for hexagonal ZnO NWs. 

4 Conclusions 

We have utilized classical molecular dynamics to study surface effects on the 
piezoelectric properties of ZnO nanowires with three different ([2110], [0110] 
and [0001]) axial orientations. A key finding is that treatment of the polar 
(0001) surface via charge compensation or surface passivation is required 
to prevent the divergence of the electrostatic energy. In the context of the 
atomistic simulations performed here, we demonstrated that not treating 
the surfaces to remove the electrostatic energy divergence results in spurious 
transformations of the initial wurtzite phase to a d-BCT phase. With regards 
to the piezoelectric properties, the piezoelectric constants of the transformed 
d-BCT phase, which occurred for nanowires with untreated surfaces, were 
nearly one order of magnitude smaller than those calculated for nanowires 
whose surfaces had been treated using either the charge compensation or 
surface passivation techniques. 

Overall, our results show that the [2110] oriented nanowires have a larger 
effective piezoelectric constant than the [0110] oriented nanowires. However, 
if proper treatment of the polar surfaces was performed, the effective piezo- 
electric constants for all nanowires were found to decrease with decreasing 
size, with all values smaller than the respective bulk ones. We further demon- 
strated the underlying atomistic mechanism for the reduction in piezoelectric 
constants, in that regardless of whether the surface expanded or contracted in 
response to surface stresses, the bond length of the Zn-0 dimer closest to the 
surface was found to decrease, thus causing a decrease in polarization at the 
nanowire surface and the corresponding reduction in effective piezoelectric 
constant. Our overall finding is therefore that due to the observed decrease 
in piezoelectric constant for all three nanowire orientations with decreasing 
size, we recommend that larger diameter square or nearly square cross section 
nanowires be utilized in practical applications if maximum energy generation 
or harvesting using ZnO nanowires is desired. 
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